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■ Recently, simple conditions for well-behaved-ness of anisotropic Voronoi diagrams have been 
proposed. While these conditions ensure well-behaved-ness of two types of practical anisotropic 
Voronoi diagrams, as well as the geodesic-distance one, in any dimension, they are both pro- 
hibitively expensive to evaluate, and not well-suited for typical problems in approximation or 
optimization. We propose simple conditions that can be efficiently evaluated, and are better 
suited to practical problems of approximation and optimization. The practical utility of this 

■ analysis is enhanced by the fact that orphan-free anisotropic Voronoi diagrams have embedded 
I triangulations as duals. 



1 Introduction and Previous Work 



■ The Voronoi diagram, or Dirichlet tesselation, is a fundamental mathematical construct, with 

i extensive practical application in diverse fields [2]. Given a discrete set of sites, its definition requires 

I a choice of distance function, which determines how the input domain is broken up into tiles of 

Q ■ points closest to a unique site. The distance typically used is the natural distance in a Euclidean 

^ . space, but generalizations based on other choices of distance can be defined [H O El [HI |12]. The 

I use of a particular distance can be motivated by the intended application. In particular, Voronoi 

(N. dia grams with respect to the geodesic distance on a Riemannian manifold can be useful in defining 

approximations of the manifold that satisfy certain properties, such as asymptotic optimality [5|ll2j. 
When considering applications in which specifying a Riemannian metric to induce a distance 
^ ■ is appropriate, the cost of computing geodesic paths for any pair of points may be very high for 

practical applications (see for instance [l3]). For this reason, fast approximations to the geodesic 
distance have been proposed [Hill], which take constant time to evaluate, provided that the metric 
can be evaluated at any point in constant time. However, while the geodesic-distance diagrams 
are guaranteed to be composed of cells, each of which contains its generating site (the diagram is 
orphan- free), this is, in general, not the case for the approximations of [6l[TT]. Since the definition 
of the dual of diagram with orphans is problematic (i.e. it may require detecting orphans and 
placing new sites inside each, a process that could potentially become recursive), and since certain 
approximation or optimization guarantees may be lost in the presence of orphans, it is natural 
to consider orphan-freedom as a basic well-behave-ness condition for these anisotropic Voronoi 
diagrams. For instance, duals of two-dimensional diagrams of the type of [6] can be shown to be 
embedded triangulations [3], if they are orphan- free. 

While the work of [6] does not provide guarantees on orphan- freedom, that of [TTj provides an 
algorithm that is guaranteed to output an orphan- free diagram, in two dimensions. More recently, 
the work of [1] provides general conditions for orphan-freedom of anisotropic Voronoi diagrams 
of both types, which hold in any number of dimensions. The condition is that the sites form an 
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(asymmetric) e-net, for a suitable value of e, which depends on the underlying metric. This condition 
is natural for problems in optimal covering, and approximation [9l [5]. However, while quite 
general, the above results have some drawbacks, which make them difficult to be applied in practice. 
In particular, for problems of average-case optimization or approximation, an optimal set of 
sites has been shown to satisfy a Delone property O [HI [TO], which is generally weaker than the 
net property. Additionally, the work of [4] uses a Lipschitz-type constant on the underlying metric 
to obtain conditions for orphan- freedom. In practice, however, computing this constant could 
potentially entail visiting all pairs of points in the domain, making it prohibitively expensive to 
compute. 

This work addresses the above drawbacks, by providing conditions for orphan-freedom of 
anisotropic Voronoi diagrams that apply to more general sets of sites, including those typically 
used for optimization and approximation, and use information from the metric that can be 
computed efficiently. 



2 Formal Setup 

Given a continuous metric Q € (a symmetric positive-definite matrix field) over a domain £ M" 
of n-dimensional Euclidean space, we begin by defining the functions 

DQ^^ia, b) = [{a - bYQbia - b)] ^''^ ,a,ben 

and 

DQ{a,b) = [{a - b)'Qa{a - b)Y^^ ,a,ben 

which are not symmetric, but are understood as distances for the purposes of constructing Voronoi 
diagrams, and where, because of their asymmetry, special attention must be given to the order of 
their arguments. 

Given a discrete set V, the Labelle-Shewchuk diagram of V [11] is composed of regions 
R{v) = {pen: Dq{p,v) < DQ{q,v)yq e ft} , v e V 



and an analogous definition of a Du/Wang diagram (6] uses the distance Dq^ instead. 

In the sequel, we will make extensive use of the square-root matrix field M, which, at each 
point, is the unique symmetric positive-definite matrix satisfying M*M = Q. As shown in [4], M 
must be of the same differentiability class as Q. Following [3], we denote the spectral matrix norm 
as /?(•), which, when applied to positive-definite matrices, returns the maximum eigenvalue, and 
define the function Pmi'), which returns the smallest eigenvalue of a positive-definite matrix. 
Finally, the (worst-case) metric variation constant (Tq is defined as 

p{M,M-^-I) 
(Jo := max ,. ^ , — -— 

a,b£n \\Ma{a-b)\\ 

where the maximum is assumed to exist (for instance, if is compact). 



3 Orphan- Free Anisotropic Voronoi Diagrams of Delone Sets 

We begin by extending the results of [Ij to sets of sites that satisfy a (C, P)-Delone property (C- 
cover, asymmetric P-packin^), which is a strict generalization of an asymmetric e-net (e-cover, 

^V^ is an asymmetric P-packing w.r.t. D if, for all v,w G V, it is either D{v, w) > P or D{w, v) > P. 
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asymmetric e packing) . The result is a simple extension of [1] , in which the proof structure is kept 
exactly as is, but some of the lemmas are slightly modified to account for the unequal relation 
between C and P. We proceed by replacing e by the cover constant C in all the technical lemmas 
of [4], except in those in which the packing property is used: Lemmas 8 and 10, which become: 



Lemma 3.1 (From Lemma 8, [1]). Let V he an asymmetric e-net w.r.t. Dq^ , and v,w G V be 
Voronoi neighbors of the resulting DW-diagram. If c G is in the Voronoi regions of v,w then 



M,{v-w)\\> P/{l + Cao) 



Lemma 3.2 (From Lemma 10, [4]). Let V be an asymmetric e-net w.r.t. Dq , and v,w ^ V be 

Voronoi neighbors of the resulting LS diagram. If k = {1 + C(Tq)/{\ — Cctq), then 



P/k < ||M^(i;-u;)|| < C{l + k) 



It is a somewhat tedious, but simple exercise to verify that, assuming a Delone property for the 
set of sites, results in the following conditions for orphan-freedom of anisotropic Voronoi diagrams. 



Theorem 3.3 (Adapted from Theorem 1, [3]). Given Q (z of worst-case variation (Jq, the 
Du/Wang diagram of a (C, P)-Delone set (with respect to Dq ) is orphan free if{P/C)'^/ (2(1 + Cao)^) 
2{Caof - ACao > 0. 



Theorem 3.4 (Adapted from Theorem 2, [4]). Given Q ^ of worst-case variation ctq, the La- 
belle/Shewchuk diagram of a {C, P)-Delone set (with respect to Dq ) is orphan free if [{P/C)'^k~'^ — 1^ — 27) /2— 
7^ - 27 > 0, where 7 = Ccjo(1 + A;) and A; = (1 + Co-o)/(l - Ca^). 

The above are similar to the analogous theorems in [4], but make explicit the relation between 
metric variation ctq, and the cover and packing constants {C,P). This generalization enables the 
application of these results not just to optimization or approximation, but also to problems in 
average-case (C^) optimization or approximation [8l[T0]. Typically, a lower metric variation means 
that a lower cover, and higher packing constants may be used while preserving orphan-freedom. If 
we fix (To, then a lower relation P/C requires a lower (stricter) cover constant to meet the well- 
behave-ness condition. For instance, for a Delone set in which P/C < 1 (less strict than a C-net, 
and similar to the one used in [lO]), in order to satisfy the conditions of the above theorems, the 
cover constant must be set lower than for the net condition {P/C = 1). Finally, note that these 
conditions reduce to the ones in [4] for the net expected. 

4 DifFerentiable and Piecewise Linear Metrics 

The main drawback of the above results, and those of [1], from a practical standpoint, is that the 
evaluation of the metric variation ctq may be prohibitively expensive in practice since, in the worst 
case, it requires visiting all pair of points in the domain. We address this problem here. 

Since ctq is a Lipschitz-type constant controlling the rate of change of Q, we may hope to find 
a simpler expression for differentiable metrics {Q £ C^) by simply taking the limit in the definition 
of (To as the pair of points becomes increasingly close, from every direction. In this way, we can 
define the differentiable-metric variation: 



(Ti := sup lim 

r^O A-i>0 



||MpAr|| 



sup lim 




= sup 



pjDrMpM-^) 
\\Mpr\\ 
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Clearly, the above expression can be evaluated by visiting each point in the domain only once. 
Note that, from the definitions of ctq and cji, it is 



o"! = sup Imi — < sup iim cjo = o"o 

In order to apply the results of Theorems l3.3l and l3.4l to the case in which we only have knowledge 
of 0"!, we can simply try to obtain an upper bound of (Tq using ui. However, while a differentiable 
function / over a compact domain in M has a Lipschitz constant L = maXpgQ [13] , there is no 
such simple relation of cJi bounding ctq from above. The reason for this is the denominator ||Mpr|| 
in the definition of ai, which prevents a simple analysis using integration and the Mean Value 
Theorem. It is, however, possible to prove similar results to those of Theorems 13.31 and 13.41 bv using 
only the differentiable-metric variation ui, without resorting to do. 

We begin by proving the following technical lemma, whose statement is the same as Lemma 2 
of [4], but replacing do by ai, but whose proof requires a new analysis using different techniques, 
and is relegated to Appendix A in the interest of conciseness. 

Lemma 4.1. If ai is the differentiable-metric variation of Q G , then for all a,b £ and every 
constant e, \\Ma{a — h)\\ < e implies 

l-eai< pmiMbM-^) < ||Mb(a - 6)||/||M,(a - < piM^M-^) < 1 + eai 

Since this is the key technical lemma from which all results in [4] stem, it can easily be ver- 
ified that the above implies well-behave-ness results using the differentiable-metric constant. In 
particular: 

Corollary 4.2. Given Q £ of differentiable-metric variation ai, the Du/Wang diagram of a 
{C,P)-Delone set (with respect to ) is orphan free if{P/Cf/ (2(1 + Caif)-2{Caif-ACai > 
0. 

Corollary 4.3. Given Q £ of differentiable-metric variation ai, the Labelle/Shewchuk diagram 
of a {C, P)-Delone set (with respect to Dq ) is orphan free if [{P/C)'^k~'^ — 1^ — 27) /2 — 7^ — 27 > 
0, where 7 = Ccji(1 + k) and A; = (1 + Ccri)/(1 - Co{). 

In particular, if the sites form an e-net, then the above results imply that a DW diagram is 
orphan-free whenever a\e < 0.09868, and a LS diagram is orphan-free if die < 0.0584. 

Relation between metric variations do and di. Interestingly, it is also possible to use 
Lemma 14.11 above to find some form of upper bound on do, using only knowledge of di. To 
this end we first need to introduce the metric variation do(C), over neighborhoods of size C as: 

do(C) := sup — — —- 

a,b£n, \\Ma{a-b)\\<C PV^a[a - 0)\\ 

which is identical to the definition of do, but only considers pairs of points, vaguely speaking, within 
range C of each other. Given Eq. [11] from the proof of Lemma |4.H the definition of di, and the 



4 



(constant velocity) parametrized straight line connecting a,b (q : [0, 1] — ?• ab), it is 



a,b£Q, 



sup 

\Ma{a- 



sup 

Ma{a- 



b)\\<C 



b)\\<C 




< 



< 



a,b£Q, 



a,beQ, 



sup 

Ma{a- 



sup 

Ma{a- 



b)\\<C 



6)||<C"'0 




(1) 



< sup O"! / 

a,ben, \Ma{a-b)\\<C JO 




<ai{l + Cai + {Cai)^/3) 



where, typically, the value of C used in the above expression will be the cover constant of the sites. 

We finally note that, given ui, it is possible to use Eq. [T]to upper bound (Tq, and plug this in 
Theorems 13.31 and 13.41 However, because it is always 1 + C(Ti + (Ccri)^/3 > 1, we can always obtain 
better (less restrictive) conditions for well-behave-ness of anisotropic Voronoi diagrams by using 
instead Corollaries 14.21 and 14.31 directly. 

4.1 Piecewise-Linear Metrics 

So far we have considered conditions for well-behave-ness of anisotropic Voronoi diagrams only for 
metrics of class and C^. In practice, however, it may be the case that the input metric is given 
as a PL function over a simplicial complex. This is typically the case in approximation problems 
in which the metric itself is derived from a previous estimate of the solution. 

In this case, the metric is almost everywhere differentiable, and therefore the integrals in Eq. [H 
and those in the proof of Lemma 14.11 can be suitably broken into pieces in which the derivative 
of M is defined. This implies that the results from Sec. [4] apply to the PL metric case, without 
change. The definition of ai, however, becomes a supremum over all points where the metric is 
differentiable, which always exists in the case of PL metrics in which M is interpolated inside 
simplicies, using values given at the vertices. In a slight abuse of notation, we replace supremums 
over the interior of simplicies with maximums, where it is understood that the supremum exists. 

In the piecewise-linear metric case, it is possible to make use of the simple structure of the 
metric to efficiently compute an upper bound of ai. Let Mf^ = d^M be the derivative of M in the 
k-th coordinate direction (which is constant inside each i-th simplex), T = {tj : i = 1, . . . ,m} the 
set of simplicies over which Q is linearly interpolated, and {vj : j € li} the set of vertices incident 
to the i-th simplex Tj. Noticing, from its definition, that cti is invariant to a scaling of r, and from 
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the convexity of p, it follows that: 



max max max 

j=l,...,m p6Ti ||r||=l 



p{DrMpMp^ 

\\Mpr\\ 



Lem.7.1 



max max max \ \\Mr,r\\pm{Mr,iDrM) 

=l,...,m pSTi ||r|| = l 



-1- 



= max max 

i=l,...,m ||r||=l 



< max max 

j=l,...,m ||r||=l 



Tmii\\Mpr\\pm{Mp{DrM) ^) 
.p&n 

mmXi{vjfp^{{DrMy^) 



Lem.7.1 p{DrM) 

= max max — — rrr 

i=l,...,m ||r|| = l mmj£/- XiiVj)'^ 

i=l,...,m ||r||=l milljg/. Ai(fj)^ 

1/2 



(2) 



- — ■ n — 

i=i,...,m mmjg/- Ai(^;j 

where Ai (p) is the smallest eigenvalue of Mp, and the numerator in the sixth line has been maximized 

by setting = p (^Mj'^ / J2i P {^'^i ) > which is clearly a unit vector. 

The bound of Eq. [2]is clearly conservative, however, crucially, it can be computed in time linear 
in the size of T, since it is a maximum over terms, one per simplex, which (for fixed dimension n) 
can be computed in constant time. 



5 Practical Considerations 

Consider the following practical scenario. A metric Q is given, which is PL over a simplicial complex 
T, in n dimensions (as in Sec. 14. 1^ by linearly interpolating the square-root M inside simplicies). 
We are also provided a finite set V of sites. 

We consider the problem of determining whether the anisotropic Voronoi diagram of V is orphan- 
free (whether of the Du/Wang, or Labelle/Shewchuk type), and if so, computing it efficiently. To 
answer this question, we can make use of Corollaries 14.21 and 14.31 To this end, we first use Eq. [2] 
to evaluate the differential-metric variation cji of Q (in time proportional to the size of T). Next, 
we may compute the packing and cover constants of V (with respect to either Dq^ or Dq ), and 
check whether the conditions of Corollaries 14.21 and 14.31 are met, in which case we are ensured that 
the anisotropic Voronoi diagram of V is orphan-free. In certain cases, such as when using the 
simple front-propagation algorithm of [3], it is the case that computing the Voronoi diagram can 
be done efficiently (in time proportional to the size of T), and the output is correct if there is a 
priori knowledge that the diagram is orphan-free. 

The missing step in the above reasoning is computing the covering and packing constants of 
V. The packing constant can be simply computed from V, and by evaluating the metric at points 
in V only, since it involves only distances between sites. Computing, or finding an upper bound 
of the cover constant, however, may require computing a full anisotropic Voronoi diagram, and 
measuring distances from all points in each Voronoi region to their generating site. In the case 
of Labelle/Shewchuk diagrams, computing the diagram, as in [11], may itself already provide an 
answer to whether the diagram is orphan- free. 
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As for Du/Wang diagrams, an efficient solution is to use the optimistic algorithm of [3]. This 
algorithm runs efficiently (in time proportional to the size of T), and produces a correct diagram 
only if the correct diagram of V is orphan-free. If it isn't, then the algorithm returns a diagram 
that differs from the true one only at the orphans. If the diagram of V isn't orphan-free, then the 
algorithm of [3] provides an upper bound of the cover constant C, which cannot satisfy Corollarv l4.21 
since the true diagram isn't orphan-free. If, on the other hand, the true diagram of V is orphan-free, 
then the algorithm provides the exact cover constant, which can be used to verify the conditions of 
Corollary 14. 2[ In conclusion, the use of the algorithm of [3], along with Corollary 14.21 to check for 
orphan-freedom, never produces a false-positive, and only produces false-negatives in cases in which 
it computes the exact cover constant, and therefore the false- negative cases can all be explained 
only by how loose the bounds of Corollary 14.21 are, and not by the fact that the algorithm of [3] is 
optimistic. We believe this point to be of practical importance. 



6 Conclusion 

Given the problem of determining whether the anisotropic Voronoi diagram of a set of sites is 
orphan-free, we have provided an analysis that results in new conditions that can be efficiently 
evaluated in the case of smooth or piecewise-linear metrics and, perhaps more importantly, which 
can be applied to typical sets of sites resulting from problems in optimization and approximation. 
While the analysis is somewhat complex and technical, the results are simple to formulate and 
verify. The emphasis throughout has been on the analysis, while the algorithmic aspects are only 
superficially discussed. We hope that this work may provide a useful piece of analysis in the 
ongoing progress in practical algorithms for the computation of anisotropic Voronoi diagrams and 
anisotropic Delaunay triangulations. 
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Appendix A 

We prove here Lemma l4.ll which states that, given a,b £ fi, and cji, such that ||Ma(a — b)\\ < e, 
then it is 

l-eai< pmiMbM-^) < ||Mb(a - 6)||/||M,(a - < p{MbM-') < I + eai 
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Clearly, it is 

\\Ma[a-b)\\ \\Ma{a-b)\\ 

and a similar argument shows that 

\\Mk{a - b)\\/\\Ma{a - b)\\ > pr,^{M^,M-^) 

We now show thats l — eai < pm{Mi,M^^), and p{Mi,M~^) < 1 + eai. We begin with the following 
simple lemma: 

Lemma 6.1. Given a non-singular matrix A € M*^^*^, it is p{A~^) = pm{A)~^ . 
Proof. If Aj, i = 1, . . . , n are the eigenvalues of A then 

p{A~^) = max |A~^| = max |Aj|~"'^ = (min |Aj|)~"'^ = pm{A)~^ 

i i i 

□ 

Consider the parametrized line segment p{\) = b + {a — b) ■ X, A € [0, 1]. Define the unit vector 
r = {a-b)/\\a - b\\, where e > ||Ma(a - 5)|| = \\a - 5|| HM^r || implies that \\a - b\\ < e/\\Mar\\. For 
any positive number k, the sub-multiplicative property of the spectral norm implies that a partition 
of [0, 1] into k equal parts satisfies: 

k-i 

i=0 

where A = 1/k. 

Taking logarithms on both sides, and the limit of A approaching 0, the above inequality becomes: 

k 

lnp(^p(A)^plo)) < |mElnp(Mp((,+i)A)M-i^)) (3) 

i=0 

l.mi:inp(/ + A^^«^^M- ) (4) 



4=0 

.^p((i+l)A) -Mp^iA) , ,_i 



lim > In 

1=0 



1 + Ap( ^^^'^-^"^ ^M-\0 



(5) 



lim V Ap( ^^((-+^)-^)-^^(-^) M7^^0 (6) 

i=0 

%iD,^aMp^,)M-')dr (7) 







rX DrMr,(r\M~\ \\M , ^r\\ 
< / p( „ ^ n^'^ ) i.i"n g'^-^ (8) 



\\Mp(r)r\\ \\Mar\ 



< 



rX 

/ ec7ip(Mp(,)M-i))dr (9) 

Equations [6] and [7] require some clarification. Equation [6] can be obtained by taking the Taylor 
expansion with remainder: ln(l + x) = x — for some constant c G [0, .x], and showing that 

the second term vanishes when taking the limit for A ^ 0. 
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Likewise, Equation [7] can be worked out backwards as: 



Since M G C^, the limit is uniformly convergent on a compact domain $7, and so we can move it 
outside of the integral. The integral exists by the continuity of its integrand, and so equals the 
limit over all partitions of [0, A] of the Riemann sum as the size of the intervals of the partition go 
to 0. Now the integrand depends on Ai, but uniformly. Therefore, for any e > there is a A' > 
such that, for all Ai < A', the Riemann sum associated with the equally-spaced partition of [0, A] 
of size Ai is within e of the true integral value. Hence 



lim /%( ^P(^A+A0-M,(.^) ^ hmf Ap( ^^«-+^^^r^^(-^^ M7^^0 



The expression 

'J 

implies that we can define an upper bound p{Mp(^x)^p(^)) ^ such that: 

/•A 

ln/x(A)= / eaifi{X)dT 
Jo 

which, taking derivatives, becomes: 

P^HX)^ = ea^p{X) (10) 

Equation [TU] is a Bernoulli ordinary differential equation. Given the initial condition //(O) = 
p{MbMj^^) = 1, it can be integrated it to obtain /i(A) = (1 — ecJiA)"^. Since /i(A) is an upper 
bound of p{Mp(^x'jM~f^Q-^), this implies 

p(Mp(;,)M-(i))<(l-eaiA)-i 
and, since a = p{l) and b = p{0), in particular 

p(M,M-i) < (l-e(Ti)-i 

By Lemma 16.11 it is 

Pm{MbM^^y' = p{MaM^^y' < (1 - eaiy' 

and thus 

Pm{MbM-^) >l-eai 

as claimed. 
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The other inequahty p{MbM^ ^) < 1 + eci can be obtained with a similar argument by using 
the identity Pm.{^ • B) > pm{A)pmiB), and the fact that Pm(^) > 1 — p{A — I): 



i=0 
k 

Km > In 



i=0 



1_ (A^Wl)A)-^p(.A) 1 



A 



A 



p(iA)> 



> 



A DrM^i^-^M"}. \\M I ^rll 
||Mp(,)r|| ^ llM^rll 

A 



where the same technical considerations as those discussed for Equations [6] and [7] apply here as 
well. 

A lower bound is //m(A) < pmiMp(^x)M~^l^) with 



In PmW 



eaiprnWdr 



which, taking derivatives, and integrating a differential equation similar to EquationllOl with initial 
condition Pm{^) = Pm{I) = !> results in Pm{^) = (1 + £criX)~^ and thus 

Pm{MaM^^) = p^{Mp^^)M-^l^) > (l + eai)-^ 

which, using Lemma l6. II implies 

p(iV4M-i)-i = p^{MaM^^) > (1 + eai)-^ 

and thus 



p{MkM~^) < 1 + ecTi 



(11) 



as claimed. 
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